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ABSTRACT 



We have compared the results of Bayesian statistical calculations and linear- 
bisector calculations for obtaining Cepheid distances and radii by the infrared 
surface brightness method. We analyzed a set of 38 Cepheids using a Bayesian 
Markov Chain Monte Carlo method that had been recently studied with a 
linear-bisector method. The distances obtained by the two techniques agree 
to 1.5% ± 0.6% with the Bayesian distances being larger. The radii agree to 
1.1% ± 0.7% with the Bayesian determinations again being larger. We interpret 
this result as demonstrating that the two methods yield the same distances and 
radii. This implies that the short distance to the LMC found in recent linear- 
bisector studies of Cepheids is not caused by deficiencies in the mathematical 
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treatment. However, the computed uncertainties in distance and radius for our 
dataset are larger in the Bayesian calculation by factors of 1.4-6.7. We give 
reasons to favor the Bayesian computations of the uncertainties. The larger un- 
certainties can have a significant impact upon interpretation of Cepheid distances 
and radii obtained from the infrared surface brightness method. 

Subject headings: Cepheids — methods: data analysis — methods: statistical 

1. Introduction 

The infrared surface brightness technique is a powerful method for determining dis- 
tances to Cepheid variables (Welch 1994, Fouque & Gieren 1997). It is independent of other 
astrophysical distance scales, nearly independent of errors in reddening, and may be applied 
to arbitrarily-chosen, individual Cepheids at great distances. However, the implementation 
of the infrared surface brightness technique and its predecessor, the visual surface bright- 
ness method, have been criticized as not mathematically rigorous in their solutions to the 
surface brightness equations (Laney & Stobie 1995, Barnes & Jefferys 1999). This leads to 
the possibility that the distances, radii, and their uncertainties may be erroneous. 

This is of current interest as Storm et al. (2004) have found a short distance to the LMC 
based on a linear-bisector analysis of six Cepheids in the LMC cluster NGC1866. They also 
found a slope to the LMC Cepheid period-luminosity relation that is substantially different 
than found in the OGLE magnitudes (Udalski et al. 1999). It is important to determine if 
their results are affected by the mathematical method used in the surface brightness analysis. 

To address the larger issue Barnes et al. (2003) developed a Bayesian Markov Chain 
Monte Carlo (MCMC) solution to the surface brightness equations that is mathematically 
rigorous. In the current paper we do a direct comparison for a significant sample of Cepheids 
between the Bayesian MCMC solution and recent linear-bisector calculations (Storm et al. 
2004, Gieren et al. 2005) to explore possible differences. The two methods used identical 
data, identical surface brightness equations and identical physical constants to ensure that 
only the mathematical approaches were compared. 

In the next section we introduce the surface brightness method for determining distances 
and radii and present the background of the infrared surface brightness method. We then 
discuss the data that are used for the two calculations. We review the linear-bisector cal- 
culations and the Bayesian MCMC calculations. In section 6 we compare the results of the 
two methods for 38 Galactic Cepheids. Finally we discuss the importance of the agreement 
and the differences that we find. 
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2. The Surface Brightness Method 

2.1. The surface brightness equations 

The infrared surface brightness method is a modification of the visual surface brightness 
technique developed by Barnes et al. (1977) and thus shares the same computational algo- 
rithm. Because solution of the equations is the important issue, we introduce the equations 
in some detail. Useful discussions of previous work have been given by Gieren, Barnes & 
Moffett (1993), Fouque & Gieren (1997), Nordgren et al. (2002), Fouque, Storm, & Gieren 
(2003), and Barnes et al. (2003). 

Barnes & Evans (1976) and Barnes, Evans, & Parsons (1976) defined a visual surface 
brightness parameter Fy as 

Fy = 4.2207 -O.IVq- 0.5 log (1) 

and also as, 

F v = \ogT e + 0.1BC, (2) 

where Vq is the stellar visual magnitude corrected for interstellar extinction, <fi is the stellar 
angular diameter expressed in milliarcseconds, T e is the effective temperature and BC is the 
bolometric correction. 

They demonstrated that Fy is well correlated with Johnson color index (V — R)q for a 
very wide range of stellar types. 

F v = A + B (V - R) o (3) 

equation (3) is called the visual surface brightness relation. 

These relations led Barnes et al. (1977) to infer a distance scale for Cepheids as follows. 
At each time t in the pulsation of the Cepheid, equations (1) and (3) may be combined to 
obtain the angular diameter variation of the star, (f>(t), 

4.2207 - 0.1V (t) - 0.5 log0(t) = A + B (V - R) (t) (4) 

In addition we infer the Cepheid's linear radius variation AR(t) about the mean radius from 
an integration of the radial velocity curve, V r (t), 

AR(t) = ~ J P ( v r(t) ~ v i) dt ( 5 ) 

where the factor — p converts observed radial velocity to the star's pulsational velocity and 
V 1 is the center-of-mass radial velocity of the star. Integration of the discrete radial velocity 
data requires that the velocity variation be appropriately modeled. 
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Substituting into equation (4), the relation among mean angular diameter </> , linear 
diameter AR(t), and distance r, we obtain 

4.2207 - 0.1V (t) - 0.5 log(0 o + 2000Ai?(t)/r) = A + B (V - R) (t). (6) 

where 0o is m milliarcseconds, r is in parsecs, and AR(t) is in AU. The factor 2000 converts 
radius to diameter and arcseconds to milliarcseconds. 

Until very recently direct solution to equation (6) for distance and diameter has proved 
daunting. Various methods have been tried, but for reasons given by Barnes & Jefferys 
(1999), none of the methods is rigorous. 

Barnes et al. (1977), Gieren et al. (1993, 1997) and Welch (1994) simplified the problem 
by solving equation (4) for 4>(t), solving equation (5) for AR(t) and then solving for r and 
0o using ordinary least-squares solution on equation (7) taking AR(t) as the independent 
variable. 

AR(t) = r(0 o + 0(t))/2OOO (7) 

However, the least-squares calculations do not properly treat the error s-in-variables problem 
that arises from uncertainty in both AR(t) and 0(t). The best that can be done using 
least-squares is the linear-bisector solution. Isobe et al. (1990) showed that linear-bisector 
performs better than other least-squares solutions when the problem is symmetric in the 
variables, as is the case here. Storm et al. (2004) adopted the linear-bisector method. 

Following Balona (1977), Laney & Stobie (1995) linearized equation (6) and applied a 
maximum-likelihood calculation to its solution. They used an iterative maximum-likelihood 
method to solve for the larger amplitude Cepheids for which the linearization is invalid. Their 
calculation addresses the errors-in-variables problem. However, results from the maximum- 
likelihood method can be quite sensitive to accurate knowledge of the uncertainties in the 
data, as discussed by Gieren et al. 1997. 

All the above methods require that the observed radial velocities be modeled in order to 
do the integration of equation (5). This creates a model selection problem. Most researchers 
choose a Fourier series to fit the radial velocities, but the number of terms to include in the 
series is subjectively chosen. 

Finally, whether least-squares or maximum-likelihood, the above calculations fail to 
treat properly the propagation of observational error through the radial velocity integral, 
equation (5). Balona (1977) introduced a widely used approximation to the error in AR(t) 
based on the uncertainty in the radial velocity data and on the assumption that the radial 
velocity data are equally spaced in pulsation phase, which is rarely true. 

It was to address these issues that the Bayesian MCMC method was developed by 
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Barnes et al. (2003). Their calculation correctly solves the errors-in-variables problem, 
objectively selects the number of terms in the Fourier series fits, and correctly propagates the 
observational error through the radial velocity integration. Moreover, the Bayesian MCMC 
calculation does not demand that a particular Fourier series fit the data; rather, each Fourier 
series that is fit to the data has a particular posterior probability. That probability is used 
as a weight in determining the other quantities sought in the solution, i.e., distance, mean 
radius, etc. 

A large number of calculations of variable star distances and radii by surface brightness 
methods are present in the literature. It is therefore of interest to determine the extent to 
which the above listed deficiencies affect those results. We will do this by comparing the 
linear-bisector calculations with Bayesian MCMC calculations. 

2.2. The infrared surface brightness technique 

Welch (1994) first showed that use of the infrared combination Kq, (V — K)q in place of 
Vq, (V — R)o in equation (6) has significant advantage in the precision of the distances and 
radii for Cepheids. He found an improvement by a factor of three in the distance uncertainty 
for the Cepheid U Sgr. Welch attributed the improved precision to several factors. First, 
the color index (V — K) is as good an indicator of surface brightness as bluer color indices 
but much less sensitive to the complications of line blanketing and surface gravity. Second, 
the Kq magnitude lies on the Rayleigh- Jeans tail of the flux distribution and is therefore less 
affected by the surface brightness variation. It is obvious that a magnitude more sensitive 
to the radius variation that is then corrected for surface brightness variation using a color 
index that is a more accurate indicator of surface brightness ought to yield superior results. 

Laney & Stobie (1995) independently examined two optical magnitude-color index com- 
binations and two infrared combinations to determine which would give the most precise radii 
of Cepheids. From both model atmosphere considerations and the precision of 49 Cepheid 
calculations, they concluded that K , (V — K) and K , (V — J)q were superior to the opti- 
cal indices. (Examination of their Tables 5-6 suggests that K , (V — K) is slightly better, 
confirming Welch's choice.) Laney & Stobie also noted that the presence of a companion 
to the Cepheid would have less effect on infrared indices than on optical ones, as Cepheids 
are more likely to have a blue main-sequence companion than a red giant companion from 
stellar evolution considerations. 

Extending Welch's (1994) work on U Sgr, Fouque & Gieren (1997) compared distance 
and radius results using V , (V — R) ; V , (V — K) ; and K , (J — K) with a somewhat 
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improved data base and improved surface brightness equations. They found the same dis- 
tances and radii from all three combinations, but much superior precision for the infrared 
color indices. The combination Vo, (V — K) seemed to be slightly preferred on the basis 
of precision. When adjusted to the same surface brightness equations, Fouque & Gieren's 
results agree within the errors with those of Welch and of Laney & Stobie. 

With a sample of 16 Galactic cluster Cepheids, Gieren, Fouque, & Gomez (1997) re- 
peated the comparison of the above three combinations. Again they concluded that the 
infrared color indices are superior to the optical in precision and also agree well with each 
other. On the other hand, the Vq, (V — R)q solutions gave distances and radii ~ 14% larger 
than the infrared solutions. While they indicated no preference between Vo, {V — K)$ and 
K , (J — K) , their Table 3 shows that V , (V — K) gave better precision than K , (J — K) 
in 13 of 16 cases. 

Based on these studies, we adopt the combination Vo, (V — K) as the basis for our 
comparison of the mathematical methods. 



Storm et al. (2004) analyzed 34 Galactic Cepheids for distances and radii using the 
Vo, (V — K)o infrared surface brightness technique. Gieren et al. (2005) enlarged the sample 
by adding four more stars. We have used the full set of 38 Cepheids. The infrared surface 
brightness relation adopted in both those studies and by us is the one determined by Fouque 
& Gieren (1997): 



The individual stellar data required for the analyses are the photometric measures 
V, (V — K), reddening E(B — V), radial velocities V r , pulsation period P, and pulsation 
phases 9. Storm et al. list the sources for photometry and radial velocities in their Table 1. 
To ensure that we used identical data in both calculations, the Bayesian analysis used the 
same input data files as used in the linear-bisector analysis. 

We added Z Lac, Y Oph, S Sge, and CS Vel to the program using data referenced in 
Gieren et al. (2005). In addition, the results for i Car given by Storm et al. were revised to 
incorporate new radial velocities, which are also referenced in Gieren et al. (2005). Again, 
the Bayesian analysis used the same input data files as used in the linear-bisector analysis 
by Gieren et al. 



3. 



The Data 



F v = 3.947 - 0.131(V-K) 
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Although the infrared surface brightness method is largely independent of the inter- 
stellar extinction, it is important that the same extinction and reddening be used in our 
two calculations. We follow the earlier studies and adopt E[B — V) from Fernie's (1990) 
tabulation for Cepheids. These E(B — V) values are listed in Table 1 along with the periods 
of the Cepheids. We adopted A v = 3.26E(B - V) and E(V — K) — 2.88E(B - V), slightly 
different than the reddening law used by Storm et al. The linear-bisector calculations were 
repeated with the new reddening law for all 38 Cepheids. 

In the integration of the radial velocity curve to obtain the linear displacements, equa- 
tion (5), the value of p is required. Some studies have adopted a single value of p for all 
Cepheids, others have adopted individual values. This choice has no effect upon our compar- 
ison of linear-bisector and Bayesian computations, provided the same value is used in both 
calculations for a specific Cepheid. Storm et al. and Gieren et al. used a relation between p 
and period developed by Gieren et al. (1989) to approximate model atmosphere results by 
Hindsley & Bell (1986): 

p = 1.39 _ O.OZlogP (9) 

where P is the pulsation period. We adopt this relation as well, using the periods given in 
Table 1 for the Bayesian MCMC calculation. 

4. The Linear-bisector Computations 

Storm et al. (2004) and Gieren et al. (2005) applied a linear-bisector, least-squares so- 
lution to the surface brightness equations. They solved equation (4) for the angular diameter 
variation (pit). The radial velocity data V r (t) were linearly interpolated in phase between the 
observed velocities and equation (5) integrated in 0.01 steps in phase to obtain the AR(t) 
variation. Finally, equation (7) was solved for r and 0o- 

To verify the solution for each Cepheid, the angular diameter variation and the displace- 
ment variation were plotted against phase in the pulsation cycle. For many of the stars, the 
angular diameter variation was a very poor match to the displacement curve in the phase 
interval 0.8 — 1.0, as illustrated in Figure 1. Storm et al. discussed the source of the poor 
fit without coming to a conclusion and decided that the best course of action was to ignore 
this phase interval in the fit. They deleted the phase interval 0.8 — 1.0 for all stars in their 
fit of equation (7) for r and 0o- 

Some Cepheids showed a small phase shift between the photometry and the radial 
velocities. The phase shift causes a loop in the upper panel and a displacement between the 
fitted curve and the points in the lower panel. For the stars showing this effect a phase shift 
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Fig. 1. — Points represent the photometrically determined angular diameters for AQ Pup. 
Crosses represent computed values that were eliminated from the fit. The line in the upper 
panel shows the bisector fit. The curve in the lower panel delineates the angular diameters 
obtained from integrating the radial velocity curve at the derived distance. 
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was determined by minimizing the scatter in the upper panel. The phase shift has already 
been imposed between the photometry and radial velocities used for Fig. 1 and thus these 
effects are not seen. The adopted phase shifts are listed in Table 1. 

There are two constants adopted in a surface brightness calculation that must be the 
same in the linear-bisector and Bayesian calculations. We chose to fix the constants in the 
Bayesian calculation to those used in the Storm et al. paper. These constants are the 
constant term in the surface brightness definition, equation (1), and the conversion factor 
from angular diameter (in milliarcseconds) and distance (parsecs) to linear radius (solar 
radii). The first of these was taken to be 4.2207 and the second as 0.10727 solar radii per 
mas-parsec (see eq. (7)). 

To perform the calculations Storm et al. used the FORTRAN subroutine SIXLIN which 
is available from Isobe et al. (1990). The computations were run on a Linux personal 
computer and took a fraction of a second per star. 

The quantities determined in the linear-bisector calculation that are of interest to us 
here are the distance r, the mean linear radius R, computed from the distance and mean 
angular diameter O , and their 1 a uncertainties. These quantities are given in Table 1. 

5. The Bayesian Markov Chain Monte Carlo Computations 

The Bayesian MCMC method that we applied here is described in detail by Barnes et 
al. (2003). It is important to realize that the Bayesian calculation treats the unknowns in 
the problem as probability distributions, not as specific values to be determined. The goal of 
the analysis is to determine the probability distribution for each parameter of interest, from 
which inferences may be drawn by appropriate means and variances. 

To compute the posterior probability distribution for each parameter requires the likeli- 
hood function on a specific model, appropriate priors, and sampling strategies for all param- 
eters. The likelihood function is the probability of obtaining the particular data given the 
model. Bayesian statistics encapsulate our understanding of the parameters in the model, 
prior to considering the data, in prior distributions. The posterior probability distribution is 
the product of the prior and the likelihood, appropriately normalized. 

The full posterior probability distribution in this problem requires solution to integrals 
in the normalization that cannot be done analytically. However, the unnormalized posterior 
probability distribution avoids these integrals and can be used to generate a Monte Carlo 
sample from the full posterior probability distribution. The techniques used to generate the 
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sample are Markov Chain techniques. 

The model for the infrared surface brightness calculation is developed by first substi- 
tuting the infrared color index for the visual color index in equation (6) and then rearrange 
as follows 

(V - K) „(t) = -J- (4.2207 - 0.1V (t) - A - 0.5 log (0 O + 2000Ai?(f)7r)) (10) 
B 

where A, B take the values given in equation (8) and where we have replaced 1/r with n, the 
parallax in arcseconds. Within this model the likelihood function is specified in a straight- 
forward way and is given in equation (12) of Barnes et al. (2003). We model the photometry 
and the radial velocity data as drawn from normal distributions with variances given by the 
observational uncertainties. Because we do not trust the quoted observational uncertainties, 
we introduce a hyper-parameter scale factor on each variance to model deviation in the scat- 
ter from that expected from the quoted uncertainties. The time variations of the photometry 
and radial velocities are modeled by Fourier series of unknown order on the pulsation phase. 

Barnes et al. discuss the priors adopted for each parameter of interest, but only one is 
relevant here. It is well-known that Cepheids are distributed within the plane of the Galaxy 
with an exponential decrease in density away from the plane. We adopted a prior on distance 
that reflects the flattened density distribution with a scale height of 70 ± 10 pc. The results 
are insensitive to reasonable changes in the scale height. 

The sampling strategy employed in this work is Markov Chain Monte Carlo using the 
Metropolis-Hastings and Gibbs algorithms, as described in Barnes et al. The art in this 
approach is to find sampling methods that explore the posterior probability distributions 
fully and efficiently Internal tests provide guidance on the completeness and efficiency of 
the sampling. All the results presented here passed those tests. In the customary manner, 
we chose a burn-in phase to improve the model selection efficiency. 

To be consistent with the linear-bisector calculation, we deleted from the Bayesian 
solution the photometry and radial velocities in pulsation phase interval 0.8 — 1.0. In the 
general model, we allow for an unknown phase shift between the photometry and the radial 
velocities; in this calculation, we fixed the phase shift at the value determined by the linear- 
bisector calculation. 

The calculations were run on a 1 GHz Macintosh G4 computer under system MacOS 
10.3 using the statistical language R-l.8.0/3 distributed by the R Development Core Team. 1 
(The R language is available for other platforms.) We used a burn-in of 1,000 samples 



See http://www.r-project.org/ . 
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followed by 10,000 samples. Because the R code is interpreted code, the calculations run 
slowly; the computations here typically took about an hour per star. 

Posterior probability distributions were determined for all the model parameters of inter- 
est. These are the parallax 7r, the mean angular diameter O , the orders of the Fourier series 
on the apparent magnitudes and the radial velocities, the aforementioned hyper-parameters 
on the observational uncertainties, the mean V^o magnitude (both intensity mean and magni- 
tude mean), and the center-of-mass radial velocity. For our purpose here, only the parallax 
and mean angular diameter are important. These were converted to distance and mean 
linear radius within the code and listed with their uncertainties in Table 1. (The erratum 
published by Barnes et al. (2003) was addressed in the computation of the radii.) 

6. Comparison of the Results 

6.1. The distances and radii 

Our goal is to determine differences in the results of the two calculations with respect 
to distance and radius. Figures 2 and 3 show the Bayesian values plotted against the linear- 
bisector values. There is no obvious difference between the values from the two calculations. 
To look for subtle differences, we computed the ratio of the Bayesian distance to the linear- 
bisector distance. A weighted least-squares fit of this ratio against log(P) gives (Figure 
4) 

Ratio(r) = 1.016(±0.007) - 0.012(±0.022)(/o#P - 1.113) (11) 

where the fit is centered on the mean period. A weighted fit for the radius ratio gives (Figure 
5) 



Ratio(R) = 1.012(±0.007) - 0.004(±0.023)(Zo#P - 1.113) (12) 

These show no evidence for any dependence on pulsation period. Therefore we take a 
weighted mean for each ratio: the distance calculations differ by 1.5% ± 0.6%, with the 
Bayesian being larger, and the radius calculations differ by 1.1% ± 0.7%, with the Bayesian 
again being larger. For comparison, the best individual distance and radius measurements 
in our dataset are for X Cyg, for which the Bayesian uncertainties are both ±2.6%, larger 
than the possible systematic difference between the Bayesian and bisector calculations. Our 
first result is that the distances and radii computed by the two methods agree quite well. 
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Bisector Distance (pc) 



Fig. 2. — The Bayesian distances are plotted against the linear-bisector distances, each with 
its 1(7 error bar. For nearby Cepheids the uncertainties are smaller than the symbols. The 
identity line is shown for reference. 
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Bisector Radius (solar radii) 

Fig. 3. — The Bayesian radii are plotted against the linear-bisector radii, each with its la 
error bar. For smaller Cepheids the uncertainties are comparable to the symbols in size. The 
identity line is shown for reference. 
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Fig. 4. — The ratio of Bayesian distance to bisector distance is plotted against logP. The 
combined uncertainty is shown as a la error bar. The line at ratio = 1.016 ± 0.006 is the 
weighted mean value. 
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Fig. 5. — The ratio of Bayesian radius to bisector radius is plotted against log P. The 
combined uncertainty is shown as a la error bar. The line at ratio = 1.012 ± 0.007 is the 
weighted mean value. 
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6.2. Uncertainties in distance and radius 

In Figures 4 and 5 the uncertainties can be seen to increase as the period is shorter. 
This is likely a result of the smaller pulsation amplitudes at shorter periods, which result in 
the photometric and velocity uncertainties having greater effect on the computed distances 
and radii. 

More importantly, a glance at Table 1, or Figures 4 and 5, shows that the uncertainties 
in the distances and radii disagree substantially between the two calculations. Typically the 
Bayesian uncertainty is more than three times the linear-bisector uncertainty. 

Because the concept of an uncertainty estimated from the Bayesian MCMC posterior 
probability distribution may not be clear, we show an example in Figure 6, the posterior 
probability distribution for the distance to U Sgr, which star is typical of our results. Over- 
plotted on the probability distribution is a normal distribution constructed for the same 
distance (592 pc), sigma (±21 pc) and area (10,000 samples). The normal distribution de- 
scribes the posterior probability distribution for the distance very well. This justifies our 
adopting the sigma of the corresponding normal distribution as a la estimator for the uncer- 
tainty in the distance (and similarly for the radius) determined in the Bayesian calculation. 
This estimator is then compared to the la estimator from the linear-bisector computation. 

We add a caveat to the previous paragraph. Because our computation determines the 
stellar parallax, not the stellar distance, (see equation (10)) the posterior probability distribu- 
tion for the distance can become asymmetric when the errors are large. As the uncertainties 
become large, the parallax posterior probability distribution becomes broad (large sigma). 
Its reciprocal, the distance posterior probability distribution, will also become broad and 
necessarily asymmetric to larger distances. The same asymmetry will arise for the radius 
posterior probability distribution in such cases because that distribution is the product of 
the angular diameter posterior probability distribution (symmetric) and the distance pos- 
terior probability distribution (asymmetric). Note that these asymmetric distributions are 
real and not mathematical artifacts; they properly represent our knowledge of the distance 
and radius, which is not true for least-squares or maximum-likelihood calculations on the 
same data. The latter methods assume symmetric errors by their very nature. Because this 
situation prevails for only a few stars in this sample, and only for stars with large errors, it 
has little effect on the weighted mean ratios of distances and radii quoted in the previous 
section. 

As we did with the distances and radii themselves, we begin by examining the behavior of 
the ratio of the Bayesian uncertainty to the linear-bisector uncertainty for the same Cepheid. 
In Figures 7 and 8 we show these ratios for the distance and radius uncertainties plotted 
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Fig. 6. — The posterior probability distribution for the distance to U Sgr is shown as a 
histogram. A normal distribution of the same mean distance, sigma and area is over-plotted. 
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against logP. Unweighted least-squares fits in these figures yield 



Ratio(r) = 3.396(±0.183) - 0.748(±0.549)(Zo#P - 1.113) 



(13) 



Ratio(R) = 3.281(±0.202) - 0.673(±0.606)(Zo#P - 1.113) 



(14) 



There is no apparent dependence of these ratios on pulsation period. Plots of the ratios 
of the distance uncertainties against distance and of the ratios of the radius uncertainties 
against radius are similarly uninformative. 

The two ratios are, however, highly correlated with each other (R = 0.99) as shown in 
Figure 9. Thus the underlying cause of the larger uncertainties in the Bayesian calculation 
is likely to be the same for the distance uncertainty and radius uncertainty. 

In section 2.1 we noted that the linear-bisector calculation does not treat the errors-in- 
variables problem rigorously nor does it properly propagate uncertainty through the radial 
velocity integration. The second of these issues will certainly lead to an underestimate of the 
uncertainties in the computed distances and radii. Because the Bayesian MCMC calculation 
does correctly address these two computational issues, we interpret the large ratio of Bayesian 
to bisector uncertainty as measuring the amount by which the linear-bisector errors have been 
underestimated. This interpretation is supported by the fact that none of the linear-bisector 
uncertainties is larger than its Bayesian counterpart. Our second result is that the linear- 
bisector calculation underestimates the uncertainties in distance and in radius substantially, 
amounting to factors of 1.4-6.7 for this dataset. This large range implies that the ratio that 
is obtained depends on the specifics of the data for the Cepheid which varies from star to 
star. 



We set out to determine whether infrared surface brightness estimates of Cepheid dis- 
tances and radii by the linear-bisector calculation are affected by the known mathematical 
shortcomings of that calculation. Based on comparison of Bayesian MCMC and linear- 
bisector calculations for 38 Cepheids using the same data, same surface brightness equa- 
tions, and same physical constants, we find that the distances and radii are not adversely 
affected but that the uncertainties in these quantities are seriously underestimated in the 
linear-bisector calculation. 

We find that Cepheid distances determined by the two calculations agree to 1.5% ±0.6% 
with the Bayesian distances being larger. This may be compared to the smallest individual 
uncertainty in distance found in the Bayesian calculation of ±2.6%. Similarly we find that 



7. 



Discussion 
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Fig. 7. — The ratio of the Bayesian distance uncertainty to the bisector distance uncertainty 
for the same Cepheid is plotted against logP. 
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Fig. 8. — The ratio of the Bayesian radius uncertainty to the bisector radius uncertainty for 
the same Cepheid is plotted against logP. 



- 21 - 




Ratio of Radius Uncertainties 



Fig. 9. — The ratio of the Bayesian distance uncertainty to the bisector distance uncer- 
tainty is plotted against the ratio of the Bayesian radius uncertainty to the bisector radius 
uncertainty. The identity line is shown for reference. 
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Cepheid radii determined by the two calculations agree to 1.1%±0.7% with the Bayesian radii 
being larger. This may be compared to the smallest individual uncertainty in radius found 
in the Bayesian calculation of ±2.6%. Any systematic difference between the mathematical 
approaches is both smaller than a la effect and smaller than the typical single-star (Bayesian) 
uncertainty dictated by the data. 

These results have an impact on interpretation of infrared surface brightness results for 
Cepheids. For example, Storm et al. (2005) used six Cepheids in an LMC cluster to infer a 
distance to the LMC of (to — M) = f8.30±0.07 mag. by means of a linear-bisector solution 
for the infrared surface brightness equations. This is less than the generally accepted distance 
of 18.50 mag. From the present work, we can say that the linear-bisector calculation used 
by Storm et al. is not the cause of the smaller distance modulus. 

Similarly, that work found a much smaller slope for the Cepheid PL relation in the 
LMC than had been found by the OGLE project. From the absence of a period dependence 
between the Bayesian MCMC distances and the linear bisector distances (Fig. 4) and from 
the overall agreement between the two, we can be certain that the smaller slope found by 
Storm et al. (2005) is not a result of using the linear-bisector, least-squares method. 

The uncertainties in distance and radius determined by the Bayesian MCMC calcula- 
tion are much larger than determined by the linear-bisector calculation. Given the known 
problems in a least-squares solution to the surface brightness equations, we interpret this 
as measuring the amount by which the linear-bisector computation underestimates the un- 
certainties. It is important to note that the ratio of Bayesian to linear-bisector uncertainty 
ranges from 1.4 to 6.7 in this set of 38 Cepheids. Clearly the amount by which the linear- 
bisector method underestimates the uncertainty depends on the specific nature of the data. 
This is expected, but inconvenient, as it is not possible to simply multiply published linear- 
bisector uncertainties by a constant correction factor. 

Our result on the underestimation of the uncertainties in a linear-bisector calculation is 
supported in the previously mentioned paper by Storm et al. (2005). As the six Cepheids 
studied by them are in an LMC cluster, we can be confident that they are at the same 
distance. The scatter in their distances is then an estimate of the true uncertainty in the 
linear-bisector distance to the cluster. This scatter was found by Storm et al. (2005) to be 
twice the formal errors of the linear-bisector distances to the Cepheids, within the range of 
results determined here. 

Other mathematical approaches have been used to solve the surface brightness equations 
for distance and radius. Ordinary least-squares assumes no error on one variable and all errors 
on the other. Inverse fit least-squares assumes the reverse. The limitations of ordinary linear 
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least-squares solutions (direct and inverse) include not only underestimation of the errors, 
but also possible systematic bias in the resulting distances and radii as discussed by Laney 
& Stobie (1995) and Gieren et al. (1997). 

The linear-bisector, least-squares calculation achieves a solution in-between the two 
other least-squares calculations (direct and inverse). Moreover, the linear-bisector error bar 
roughly corresponds to the difference between the two least-squares solutions. Maximum 
likelihood uses information on errors on one or both variables to choose a result between the 
two results of linear least-squares (direct and inverse fits). As a result, maximum likelihood 
cannot differ by more than one linear-bisector sigma from ordinary least-squares and even 
less from a linear-bisector fit. As we have shown that this linear-bisector sigma is about 
1/3 of the Bayesian sigma, the maximum likelihood results should also be about 1/3 of the 
Bayesian values. Our results for the linear-bisector solutions thus suggest that the maximum 
likelihood method would yield distances and radii in the infrared surface brightness method 
that are unbiased if the uncertainties in the data are well understood. Unfortunately, these 
uncertainties are often not well understood. First, Barnes et al. (2003) showed that quoted 
uncertainties in Cepheid photometry and radial velocities are usually underestimated. Sec- 
ond, maximum-likelihood calculations usually adopt an approximation for the uncertainty 
in the displacements advocated by Balona (1977) for equally spaced velocity data. Our work 
demonstrates that this approximation does not apply to typical unequally-spaced radial ve- 
locity curves; if it did apply, the linear-bisector method would have yielded uncertainties in 
distance and radius close to those of the Bayesian MCMC calculation. Thus we expect the 
maximum likelihood uncertainties to be underestimated as are the linear-bisector uncertain- 
ties and for the same reasons. 
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Table 1. Distance and Radius Results 



Cepheid 


Period 


E(B-V) 


A9 


r 




r 




R 




R 






days 


mag 
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pc 


pc 


pc 


pc 




RlT) 




RlT) 




LSQ 




Bayes 




LSQ 




Bayes 




SU Cas 


1.949322 


0.287 


0.000 
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±14 


460 


±36 


29.4 


±0.9 


31.8 


±2.4 


EV Set 


3.090990 


0.679 


0.045 


1714 


±83 
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±453 


34.2 


±1.6 


36.1 


±9.0 


BF Oph 
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0.035 
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±11 
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±34 


32.0 


±0.5 
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±1.4 


T Vel 
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0.281 


0.000 
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±25 
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33.6 


±0.9 


40.0 
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5 Cep 


5.366341 
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±5 
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±20 


42.0 


±0.9 


43.1 


±3.2 


CV Mon 


5.378793 


0.714 


0.015 


1527 


±24 


1677 


±146 


40.6 


±0.6 


44.6 


±4.0 


V Cen 


5.49392 


0.289 


0.000 
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±20 


636 


±48 


42.0 


±1.2 


39.6 


±2.9 


CS Vel 


5.904740 


0.847 


-0.005 


2908 


±86 


3204 


±297 


38.5 


±1.1 


42.3 


±4.0 


BB Sgr 


6.63699 


0.284 


-0.035 


789 


±10 


729 


±43 


49.9 


±0.6 


45.6 


±2.6 


U Sgr 
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0.000 
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±6 
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47.7 


±0.5 


49.1 
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±5 
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48.4 


±1.0 


49.2 


±2.8 


SSge 


8.382086 


0.127 


-0.010 
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±11 


713 


±25 


59.1 


±1.0 


60.7 


±2.2 


S Nor 


9.754244 


0.189 


0.000 


949 


±14 
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±51 


70.7 


±1.0 


73.8 


±3.7 


Z Lac 
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0.404 


-0.005 


1905 


±38 
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±55 


71.4 


±1.4 


73.7 


±2.1 
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10.95337 


0.260 


-0.040 
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±17 


1601 


±90 


69.4 


±0.7 


67.2 


±3.7 
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11.287 


0.315 


0.000 
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±142 
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±420 


67.1 


±5.7 


79.2 


±16.9 
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11.63641 


0.413 


-0.005 


3224 


±125 


3227 


±433 


74.1 


±2.9 


73.5 


±9.9 


U Nor 


12.643710 


0.892 


0.000 


1325 


±37 


1310 


±68 


76.3 


±2.1 


75.0 


±3.9 


BN Pup 


13.67310 


0.438 


0.000 


3802 


±88 
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±277 


83.3 


±1.9 


90.9 
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14.146400 


0.478 


0.000 
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±129 
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±438 


90.2 


±2.3 


90.3 


±7.8 
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15.03618 


0.448 


0.000 


3550 


±64 
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±186 


86.7 


±1.6 


91.2 


±4.5 



Table 1 — Continued 
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X Cyg 


16.386332 


0.288 


0.000 


1194 


±10 


1233 


±33 


105.5 


±0.9 


108.6 


±2.8 


Y Oph 


17.12413 


0.655 


-0.035 


573 


±8 


590 


±42 


89.5 


±1.2 


92.1 


±6.6 


VY Car 


18.9154892 


0.243 


-0.020 
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±20 
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±98 


112.7 


±1.1 


112.2 


±5.6 


RY Sco 


20.320144 


0.777 


0.000 


1215 


±19 


1253 


±62 


99.9 


±1.6 


103.1 


±5.1 


RZ Vel 


20.3969 


0.335 


0.000 
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±21 
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±64 


114.8 


±1.5 
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±4.6 
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0.467 
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±38 
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±139 
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±2.6 
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23.01320 
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0.000 
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±114 
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±177 
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±3.4 


108.5 


±5.2 


VZ Pup 


23.17100 


0.471 


0.000 


4027 


±104 
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±243 


96.9 


±2.5 


95.5 


±5.9 


SW Vel 


23.443130 


0.349 


-0.020 
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±28 


2557 


±85 


117.3 


±1.4 


122.2 


±4.1 
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27.03428 


0.209 


0.000 
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±36 
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±55 


149.6 


±3.8 
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±5.7 


RY Vel 


28.129259 


0.562 


-0.005 


2458 


±36 


2410 


±185 


139.9 


±2.1 


136.4 


±10.5 


AQ Pup 


30.1040 


0.512 


-0.055 


3106 


±64 


3251 


±216 


147.8 


±3.1 


154.8 


±10.4 


KN Cen 


34.029641 


0.926 


0.005 


4011 


±83 


3990 


±231 


185.8 


±3.9 


184.3 


±10.8 


i Car 


35.54804 


0.170 


-0.025 


561 


±6 
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±19 


179.0 


±1.8 


179.9 


±6.4 


U Car 


38.81234 


0.283 


-0.035 


1497 


±17 


1545 


±68 


157.6 


±1.9 


161.8 


±7.2 


RS Pup 


41.4400 


0.446 


0.000 
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±59 


1984 


±122 


208.1 


±6.1 


205.1 


±12.8 


SV Vul 


44.994772 


0.570 


-0.045 


2553 


±43 


2438 


±149 


224.0 


±3.8 


213.7 


±13.8 
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